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Solar activity is studied using a cluster analysis of the time-fluctuations of the sunspot number. 
It is shown that in an Historic period the high activity components of the solar cycles exhibit 
strong clustering, whereas in a Modern period (last seven solar cycles: 1933-2007) they exhibit 
a white-noise (non-) clustering behavior. Using this observation it is shown that in the Historic 
period, emergence of the sunspots in the solar photosphere was strongly dominated by turbulent 
photospheric convection. In the Modern period, this domination was broken by a new more active 
dynamics of the inner layers of the convection zone. Then, it is shown that the dramatic change of 
the sun dynamics at the transitional period (between the Historic and Modern periods, solar cycle 
1933-1944yy) had a clear detectable impact on Earth climate. A scenario of a chain of transitions 
in the solar convective zone is suggested in order to explain the observations, and a forecast for the 
global warming is suggested on the basis of this scenario. A relation between the recent transitions 
and solar long-period chaotic dynamics has been found. Contribution of the galactic turbulence 
(due to galactic cosmic rays) has been discussed. These results are also considered in a content of 
chaotic climate dynamics at millennial timescales. 

PACS numbers: 92.70.Qr, 92.70.Mn, 96.60.qd, 98.70.Sa 



INTRODUCTION 

The sunspot number is the main direct and reliable 
source of information about the solar dynamics for his- 
toric period. This information is crucial, for instance, for 
analysis of a possible connection between the sun activity 
and the global warming. In a recent papers [l|,[2] results 




FIG. 1: Sunspot number (SSN, monthly) vs time [H]. The 
dashed straight lines separate between periods of different in- 
tensity of the solar activity. The solid curve shows the global 
temperature anomaly (combined land and marine, 7-years 
running average) [20l ]. 



of a reconstruction of the sunspot number were presented 
for the past 11,400 years. The reconstruction shows that 
01 : " . .the level of solar activity during the past 70 years 
is exceptional, and the previous period of equally high 
activity occurred more than 8,000 years ago" (section III 
in Fig. 1) . 

Figure 1 gives first indication ('by eye'), that the 
prominent maxima of the global temperature data (solid 
curve in the Fig. 1) correspond to transitions between pe- 
riods of different intensity of the Sun activity (character- 
ized by the monthly sunspot number-SSN) . This observa- 
tion can be considered as an indication of a strong impact 
of the solar activity transitions on the global Earth cli- 
mate. Therefore, understanding of the physical processes 
in the Sun, which cause these activity transitions, seems 
to be crucial for any serious forecast for global Earth cli- 
mate. 

SUNSPOTS 

When magnetic field lines are twisted and poke 
through the solar photosphere the sunspots appear as 
the visible counterparts of magnetic fiux tubes in the 
convective zone of the sun. Since a strong magnetic field 
is considered as a primary phenomenon that controls gen- 
eration of the sunspots the crucial question is: Where has 
the magnetic field itself been generated? The location of 
the solar dynamos is the subject of vigorous discussions 
in recent years. A general consensus had been developed 
to consider the shear layer at the bottom of the convec- 
tion zone as the main source of the solar magnetic field 
[1] (see, for a recent review [1]). In recent years, however, 
the existence of a prominent radial shear layer near the 
top of the convection zone has become rather obvious and 
the problem again became actual. The presence of large- 
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scale meandering flow fields (like jet streams), banded 
zonal flows and evolving meridional circulations together 
with intensive multiscale turbulence shows that the near 
surface layer is a very complex system, which can sig- 
nificantly affect the processes of the magnetic field and 
the sunspots generation. There could be two sources for 
the poloidal magnetic field: one near the bottom of the 
convection zone (or just below it Q), another resulting 
from an active-region tilt near the surface of the convec- 
tion zone. For the recently renewed Babcock-Leighton 
solar dynamo scenario, for instance, a combina- 
tion of the sources was assumed for predicting future so- 
lar activity levels 0, [I]- In this scenario the surface 
generated poloidal magnetic field is carried to the bot- 
tom of the convection zone by turbulent diffusion or by 
the meridional circulation. The toroidal magnetic field is 
produced from this poloidal field by differential rotation 
in the bottom shear layer. Destabilization and emer- 
gence of the toroidal fields (in the form of curved tubes) 
due to magnetic buoyancy can be considered as a source 
of pairs of sunspots of opposite polarity. The turbulent 
convection in the convection zone and, especially, in the 
near-surface layer captures the magnetic fiux tubes and 
either disperses or pulls them trough the surface to be- 
come sunspots. 

The magnetic field plays a passive role in the photo- 
sphere and does not participate significantly in the tur- 
bulent photospheric energy transfer. On the other hand, 
the very complex and turbulent near-surface layer (in- 
cluding photosphere) can significantly affect the process 
of emergence of sunspots. The similarity of light elements 
properties in the spot umbra and granulation is one of the 
indication of such phenomena. 

The present paper reports a direct relation between 
the fluctuations in sunspot number and the temperature 
of photospheric turbulent convection in the Historic pe- 
riod. This relation allows for certain conclusions about 
the generation mechanisms of the magnetic fields and the 
sunspots. In the Modern period (section III in Fig. 1, 
cf. [Il,[l]) the relative role of the surface layer (photo- 
sphere) in the process of emergence of sunspots decreased 
in comparison with the Historic period (section II in Fig. 
1), implying a drastic increase of the relative role of the 
inner layers of the convection zone in the Modern period. 



TIME-CLUSTERING OF FLUCTUATIONS 

In order to extract new information from sunspot num- 
ber data we apply the fluctuation clustering analysis sug- 
gested in the Ref. 



10| (see also Ref. [ll|). For a time 
depending signal we count the number of 'zero'-crossings 
of the signal (the points on the time axis where the signal 
is equal to zero) in a time interval r and consider their 
running density Ut-- Let us denote fluctuations of the 
running density as Sn-r — Ur — {nT), where the brackets 
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FIG. 2: The standard deviation for Stlt vs r (Historic period) 
in log-log scales. The straight line (the best fit) indicates the 
scaling law Eq. (1). 
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FIG. 3: The standard deviation for Srir vs r (Modern period) 
in log-log scales. The straight line (the best fit) indicates the 
scaling law Eq. (1). 



mean the average over long times. We are interested in 
scaling variation of the standard deviation of the running 
density fluctuations (Jn^)^/^ with t 
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For white noise signal it can be derived analytically 
that a = 1/2 (see also [l3|). The same con- 
sideration can be applied not only to the 'zero'-crossing 
points but also to any level-crossing points of the signal. 

One can see that for the Modern period (section HI in 
Fig. 1) the solar activity is significantly different from 
that for the Historic period (section II in Fig. 1). There- 
fore, in order to calculate the cluster exponent (if exists) 
for this signal one should make this calculation separately 
for the Modern and for the Historic periods. We are in- 
terested in the active parts of the solar cycles. Therefore, 
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FIG. 4: The standard deviation for Sut vs r for the tempera- 
ture fluctuations in the Rayleigh-Bernard convection labora- 



tory experiment with Ra = 1.4 x 
indicates the scaling law Eq. (1). 
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FIG. 5: The standard deviation for StIt vs. r for radial com- 
ponent Br of the interplanetary magnetic field, as measured 
by the ACE magnetometers for the last solar cycle (hourly 
average I26|]). 



for the Historic period let us start from the level SSN=85. 
The set of the level-crossing points has a few large voids 
corresponding to the weak activity periods. In the tele- 
graph signal Ref. [10], corresponding to the data set, 
the large voids became so clear detectable that there is 
no problem to cut them off. Then, the remaining data 
have been merged providing a statistically stationary set 
(about 10"* data points). The robustness of the procedure 
has been successfully checked. 

Fig. 2 shows (in the log-log scales) dependence of the 
standard deviation of the running density fluctuations 
{Sn'^y^'^ on r for this data set. The straight line is drawn 
in this figure to indicate the scaling (1). The slope of 
this straight line provides us with the cluster-exponent 
a = 0.37 ± 0.02. This value turned out to be insensitive 
to a reasonable variation of the SSN level. Results of 
analogous calculations (including the large voids 'cut off' 
procedure) performed for the Modern period are shown 
in Fig. 3 for the SSN level SSN=125. The calculations 
performed for the Modern period provide us with the 
cluster-exponent a = 0.5 ± 0.02 (and again this value 
turned out to be insensitive to a reasonable variation of 
the SSN level). 

The exponent a ~ 0.5 (for the Modern period) indi- 
cates a random (white noise like) situation. While the 
exponent a ~ 0.37 (for the Historic period) indicates 
strong clustering. The question is: Where is this strong 
clustering coming from? It is shown in the paper [10| 
that turbulence produces signals with strong clustering. 
Moreover, the cluster exponents for these signals depend 
on the turbulence intensity and they are nonsensitive to 
the types of the boundary conditions. Fortunately, we 
have direct estimates of the value of the main parame- 
ter characterizing intensity of the turbulent convection 
in photosphere: Rayleigh number Ra ~ 10^^ (see, for in- 



stance [15|). In Fig. 4 we show calculation of the cluster 
exponent for the temperature fluctuations in the classic 
Rayleigh-Bernard convection laboratory experiment for 
Ra ~ 10"'^"'^ (for a description of the experiment details 
see [iBl). The calculated value of the cluster exponent 
a = 0.37 ± 0.01 coincides with the value of the cluster- 
exponent obtained above for the sunspot number fluctu- 
ations for the Historic period. The value of the Rayleigh 
number Ra in the photosphere for the Historic period 
has the same order as for the Modern period: Ra 10^^ 
(see next Section). Therefore, the photospheric temper- 
ature fluctuations can produce the strong clustering of 
the sunspot number fluctuations for the Historic period. 
This seems to be natural for the case when the photo- 
spheric convection determines the sunspot emergence in 
the photosphere. However, in the case when the effect of 
the photospheric convection on the SSN fluctuations is 
comparable with the effects of the inner convection zone 
layers on the SSN fluctuations the clustering should be 
randomized by the mixing of the sources, and the cluster 
exponent a ~ 0.5 (similar to the white noise signal). The 
last case apparently takes place for the modern period. 

Since the Rayleigh number Ra of the photospheric con- 
vection preserves its order Ra 10^^ with transition 
from the Historic period to the Modern one (see next 
Section), we can assume that just significant changes of 
the dynamics of the inner layers of the convection zone 
(most probably - of the bottom layer) were the main rea- 
sons for the transition from the Historic to the Modern 
period. 
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FIG. 6: Cumulative aa-index vs. time. Daily aa-index was FIG. 7: The standard deviation for Sut vs. t (in log-log 
taken from [3(|]. The arrow indicates beginning of the transi- scales) for the daily aa-index ^3]: circles correspond to the 
tional solar cycle. Historic period and crosses correspond to the Modern period. 



MAGNETIC FIELDS IN THE SOLAR WIND 
AND ON EARTH 

Although for the Modern period the turbulent con- 
vection in the photosphere has no decisive impact on 
the sunspots emergence, the large-scale properties of the 
magnetic field coming through the sunspots into the pho- 
tosphere and then to the interplanetary space (so-called 
solar wind) can be strongly affected by the photospheric 
motion. In order to be detected the characteristic scal- 
ing scales of this impact should be larger than the scal- 
ing scales of the interplanetary turbulence (cf. [23|.[25t). 
In particular, one can expect that the cluster-exponent 
of the large-scale interplanetary magnetic field (if ex- 
ists) should be close to a 0.37. In Fig. 5 we show 
cluster-exponent of the large-scale fluctuations of the ra- 
dial component Br of the interplanetary magnetic field. 
For computing this exponent, we have used the hourly 
averaged data obtained from Advanced Composition Ex- 
plorer (ACE) satellite magnetometers for the last so- 
lar cycle [2^. As it was expected the cluster-exponent 
a ~ 0.37 ±0.02. Analogous result was obtained for other 
components of the interplanetary magnetic field as well. 

At the Earth itself the solar wind induced activity is 
measured by geomagnetic indexes such index (in 

units of 1 nT). This, the most widely used long-term ge- 
omagnetic index [S^ , presents long-term geomagnetic ac- 
tivity and it is produced using two observatories at nearly 
antipodal positions on the Earth's surface. The index is 
computed from the weighted average of the amplitude 
of the field variations at the two sites. It was recently 
discovered that variations of this index are strongly cor- 
related with the global temperature anomalies (see, for 
instance, [2^,[2§]). In this paper, however, we will be 
mostly interested in the clustering properties of the aa- 



index and their relation to the modulation produced by 
the photospheric convective motion. The point is that 
the data for the aa-index are available for both the Mod- 
ern and the Historic periods. This allows us to check the 
suggestion that the Rayleigh number of the photospheric 
convection has the same order for the both mentioned 
periods. The transition between the two periods can be 
seen in Figure 6, where we show the cumulated aa-index 
ait) 

Saa{t)= I a[t')dt' (2) 
Jo 

The arrow in this figure indicates beginning of the tran- 
sitional solar cycle. 

Figure 7 shows cluster-exponents for both the Modern 
and the Historic periods calculated for the low intensity 
levels of the aa-index (for the Historic period the level 
used in the calculations is aa-index=10nT, whereas for 
the Modern period the level is aa-index=25nT, the aa- 
index was taken from [lO]). The low levels of intensity 
were chosen in order to avoid effect of the extreme phe- 
nomena (magnetic storms and etc.). The straight line in 
Fig. 7 is drawn to show the expected value of the cluster- 
exponent a ~ 0.37 for the both periods. It means that 
indeed for both the Historic and the Modern periods the 
Rayleigh number Ra ~ 10^^ (cf. Fig. 4) in the solar 
photosphere. 

IMPACT OF THE SOLAR DYNAMICS ON 
EARTH CLIMATE 

The transitions of the solar convection zone dynamics 
can affect the earth climate through (at least) two 
channels. First channel is a direct change in the heat 
and light output of the Sun, especially during the 
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FIG. 8; Global temperature anomaly (solid curve, combined 
land and marine) ^C?| and atmospheric CO2 (circles) '2?| vs 
time. Relevant daily SSN are also shown. 
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FIG. 9: Comparison of the CO2 growth for the periods 1900- 
1930yy [H and 1977-2008yy [H]. The soUd straight lines 
(best fit) indicate a linear growth. 
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FIG. 10: Comparison of the temperature anomaly (combined 
land and marine) ^20;) for the periods 1912-1942yy and 1977- 
2007yy. 



8 shows global temperature anomaly (the solid curve) 
2C| and atmospheric CO2 mixing ratios (the circles) 
2l| vs. time. Relevant SSNs are also shown. The 
dashed straight line (corresponding to the peak B) 
separates between Historic and Modern periods. One 
can see, that the transitional solar cycle (1933-1944yy) 
is characterized by dramatic changes both in the global 
temperature (a huge peak) and in the atmospheric 
CO2 (complete suppression of the growth). After 
the transitional period one can observe (during the 
three solar cycles: 1944-1977yy) certain growth in the 
temperature anomaly and an unusually fast growth in 
the atmospheric C02- These three solar cycles seems 
to be an aftershock adaptation of the global climate to 
the new conditions. Then, starting from 1977 year, the 
atmospheric CO2 growth returns to its linear trend as 
before the transition B, but now the rate of the growth 
is about four times larger than before the transition (sec 
Fig. 9). On the other hand, it can be seen from Fig. 1 
that the temperature anomaly before the maximum C 
returns to about the same growth pattern as it was just 
before the transition B. A quantitative comparison of 
these patterns, shown in Figure 10, indicates only about 
20% difference in the growth rate (cf. also Refs. [3l|, 
[s^p. The thirty year periods used for this comparison 
could find a support in the Appendix B. 



transitional solar cycle 1933-1944yy. Second channel 
is related to the strong increase of the magnetic field 
output in the interplanetary space through the sunspots. 
The interplanetary magnetic field interacts with the 
cosmic rays. Therefore, the change in the magnetic 
field intensity can affect the Earth climate through the 
change of the cosmic rays intensity and composition 
(see, for instance, Refs. 17|.[l8|.[l9t) 



more closely to the transition B in Fig. 



Let us look 
1. The Fig. 



If one consider the impact mechanism, related to the 
energetic cosmic (charged) particles and the correspond- 
ing screen effect of the magnetic fields, one could expect 
that this impact mechanism should be also sensitive to 
the 11-year solar cycle variability (though considerably 
less than to the transitional effects). To find fingerprints 
of this variability in the global temperature anomaly data 
is not a trivial task due to comparatively (to the 11- 
year cycle) short period of observations and due to the 
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FIG. 11; The solid curve corresponds to energy spectrum 
of the fluctuating (noise-hke) component of the decomposed 
global temperature anomaly time series (combined land and 
marine 20]). The dashed curve corresponds to energy spec- 
trum of fluctuating component of the SSN data. 



FIG. 12: Energy spectrum of the AR (fluctuating) component 
of the decomposed global temperature anomaly time series 
(combined land and marine [20|]). 



statistically non-stationary character of these data (with 
considerable time trends). To extract these fingerprints 
we will use a time series decomposition method devel- 
oped in [1^ , [131 , [i^l ■ A time series is decomposed into a 
trend (like that shown in Fig. 1 as the solid curve) and 
a fluctuating (noise-like) components (see Appendix A). 
The method applies state space modeling and Kalman 
filter. Parameters are estimated by maximum likelihood 
method. Obtained by this way the fiuctuating (noise- 
like) component of the (filtered) signal presents a sta- 
tistically stationary time series and, therefore, it can be 
subject of a spectral analysis. Result of this analysis 
for the glob al temperature anomaly (combined land and 
marine [20|) is shown in Figure 11. Although we have 
deal with a comparatively short data series one can clear 
see the main spectral peak corresponding to the (about) 
11-year solar cycle impact (cf. Refs. [3| ) IIS j ISl ) • 
one also introduces an autoregressive (AR) component 
into the time series decomposition (see Appendix A) to 
describe the fluctuating part of the global temperature 
anomaly, then one can compute energy spectrum of the 
AR component. Such spectrum is shown in figure 12. 
One can clear see the spectral peak corresponding to the 
solar 11-year cycle impact. 

It is also interesting to look at the reduction in arc- 
tic ice extent. Figure 13 shows a trend component of 
the anomalies of ice August extent in the Chukchi sea 
(the data are taken from [s^) obtained by the time se- 
ries decomposition. The dashed straight line indicates 
beginning of the transitional solar cycle. Insert to Fig. 
13 shows spectrum of corresponding fiuctuating (noise- 
like) component of the (filtered) signal (cf. Fig. 11). 
The dashed straight line indicates beginning of the tran- 
sitional solar cycle. This arctic sea is located sufficiently 
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FIG. 13: Trend of anomalies of ice August extent in the 
Chukchi sea [361 ]. The insert shows spectrum of fluctuating 
(noise-like) component of the time series decomposition. 

far from the North Atlantic to diminish its influence, 
whereas influence of the Pacific ocean on the Chukchi 
sea local climate is not very significant. 

A FORECAST 

At present time we can already observe a peak C in the 
global temperature anomaly (see Fig. 1). The last obser- 
vation can be considered as an indication of the current 
transition to a new section IV in the solar activity. The 
fact that in the 'last' period (1977-2007yy) the temper- 
ature anomaly growth returned to the same pattern as 
just before the transition peak B (see Fig. 10) provides 
an additional support to this suggestion. 
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FIG. 14: The upper solid curve corresponds to the observed 
yea rly global temperature anomaly data (7-year smoothed) 
[20I ]. The lower solid curve corresponds to the observed SSN 
data (year-smoothed) [ij ■ The dashed curves are an ARMA- 
prediction: for SSN and for global temperature anomaly. 
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FIG. 16: Spectrum of the galactic cosmic ray count rate fluc- 
tuations in the in-ln scales (the reconstructed data for pe- 
riod 1611-2007yy have been taken from Q)- The straight 
dashed line is drawn to indicate the Kolmogorov scaling law 
i5(/)~/(-5/3). 



1 1 - year cycle 




-6 -5 -4 -3 -2 -1 



In / [y-1] 

FIG. 15: ARMA spectrum of the global temperature anomaly 
fluctuations in the in-ln scales. The straight dashed line is 
drawn to indicate the Kolmogorov scaling law E{f) ^ j(~5/3) 



To make a forecast for the new section IV let us as- 
sume that there was a storage of certain part of the en- 
ergy coming from the radiative zone during the period of 
the section II (Fig. 1). Then, in the period of the section 
III this stored energy was released and transported to 
the Sun surface. This release was not in the form of the 
thermal convection, but in the form of strongly concen- 
trated toroidal magnetic field structures, which results in 
the intense emergence of the sunspots in the active parts 
of the Modern period (section III in Fig. 1). It is inter- 
esting that the storage and release periods have about 
the same length: ~ 60 years (see also Appendix B). The 



storage of energy in a growing magnetic field at the bot- 
tom of the convection zone or in the overshoot (by the 
convective downflows) area just under the bottom can be 
considered as the most plausible for this scenario. Con- 
vective overshooting provides means to store magnetic 
energy below the convection zone by resisting the buoy- 
ancy effect for large time scales. Level of subadiabaticity 
in the overshoot layer, needed in order to store the mag- 
netic flux for sufficiently long times, can be significant. 
Expected intensity of 10^ G results an energy density 
that is in order of magnitude larger than equipartition 
value. Though, suppression of convective motion by the 
magnetic flux tubes themselves can be considered as an 
additional factor allowing to achieve the sufficient level 
of subadiabaticity (see, for instance . [3| . [23l| ) . 

Since the section III in the solar activity was a 
release one, it is naturally to assume that the next 
section (IV) in the solar activity will be a 'storage' 
one, similar to the 'storage' section II. As for the 
corresponding Earth temperature anomaly, the section 
IV is also expected to be similar to the section II. But 
the temperature anomaly from the beginning of the 
section IV will be shifted upward on about 0.6°C in 
comparison with the section II. It then will be not 
surprising if the total length of the section IV also will 
be about 60 years (cf. Appendix B) and the section 
will be finished with a maximum D higher than 0.6°C. 
A very schematic presentation for a first part of the 
section IV is shown in figure 14 (see Appendix A for an 
Auto-Regressive-Moving- Average (ARMA) model used 
for the prediction, cf. also Ref. [39|). Figure 15 shows 
corresponding ARMA model spectrum of the global 
temperature anomaly fluctuations (in In-ln scales). A 
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FIG. 17; Spectrum of the global temperature anomaly fluc- 
tuations in the In-ln scales for a cutted observed base ARMA 
model. The straight dashed line is drawn to indicate the Kol- 
mogorov scaling law E{f) ^ j(~5/3) 



improving the paper. A software provided by K. Yosh- 
ioka was used at the model computations. 



APPENDIX A 

The above applied time series decomposition method 
was developed in the Refs. [11] , 34 1 , [sl] . In this method a 



nonstationary mean time series is decomposed into trend, 
autoregressive and noise components. A random process 
Xt in discrete time t is defined by the following expres- 
sion: 



Xt = Tt+et + AR 



(m) 



(Al) 



where tj is a trend component, St is an observational 
noise component, and 



(m) _ 



ARr> = 



(1) 



{A2) 



dashed straight line in this figure indicates a scaling 
with '-5/3' exponent: E{f) ~ /"^/^ Although, the 
scaling interval is short, the value of the exponent is 
rather intriguing. This exponent is well known in the 
theory of fluid (plasma) turbulence and corresponds to 
so-called Kolmogorov's cascade process (see, for instance 
[4o|). This process is very universal for turbulent fluids 
and plasmas 41|,[43|. For turbulent processes in Earth 
and in Heliosphere the Kolmogorov-like spectra with 
such large time scales cannot exist. Therefore, one 
should think about a Galactic origin of Kolmogorov 
turbulence (or turbulence-like processes [iH) with such 
large time scales. This is not surprising if we recall 
possible role of the galactic cosmic rays for Earth climate 
(see, for instance, 17|, 18|,[l^). In order to support 
this point we show in figure 16 spectrum of galactic 
cosmic ray intensity at the Earth's orbit (reconstruction 
for period 1611-2007yy 0, cf. also Ref. ^). One 
can compare Fig. 16 with Figs. 15 corresponding to 
the global temperature anomaly fluctuations (see also 
Appendix C). Presence of the turbulent-like component 
in the low-frequency fluctuations makes any long-range 
prediction a very difficult task (see also Appendix B on 
a chaotic element in solar dynamics). A chaotic element 
in global climate dynamics itself (Appendix C) presents 
another problem for the forecast. Such internal forcing 
agents as volcanoes and anthropogenic greenhouse gases, 
for instance, have a potential to change drastically the 
future dynamics of the global climate. 

The author is grateful to J.J. Niemela, to K.R. Sreeni- 
vasan, and to SIDC-team, World Data Center for the 
Sunspot Index, Royal Observatory of Belgium for shar- 
ing their data and discussions. The author is also grateful 
to C.H. Gibson and to R.A. Treumann for their help in 



is a globally stationary autoregressive (AR) component 
of order m (a^ are the coefficients of a recursive filter, 
and ry^^-* is a Gaussian white noise). 

A fc-th order stochastically perturbed difference equa- 
tion describes evolution of the trend component 



Tt 



,(2) 



Vt 



(A3) 



where ijl is a Gaussian white noise. 

In the Ref. algorithms, using state space modeling 
and Kalman filter, are suggested for estimation of the 
model parameters by the maximum likelihood method. 

Using so-called autoregressive moving average 
(ARMA) model one could try to make a prediction. For 
this model a random process Xf in discrete time t is 
defined by the following expression: 



(A4) 



where et is Gaussian white noise with zero mean 
(3^ , 34 1 , [3^ . Again, state space modeling and Kalman 



filter, are used for estimation of the model parameters 
by the maximum likelihood method. Figure 14 shows 
an optimistic prediction produced with this model for 
a first part of the section IV. Though, continuation of 
the analogy to the section II makes the scenario much 
less optimistic. A quasi-Newton method and dynamic 
restrictions were used for optimization of the model pa- 
rameters. 

We can also play with the model by cutting the ob- 
served time series. This procedure reduces relative con- 
tribution of the low-frequency fluctuations, presumably 
of Galactic origin (cf. Fig. 17 and Fig. 15). Then, the 
insert to the Fig. 17 shows less steep decline of the global 
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temperature anomaly at this (cutted observed base) ver- 
sion of the ARMA model prediction. This game has also 
another end. One can assume that using a more long 
time series (observed base), than that we have used in 
computing the results shown in Fig. 14, one could ob- 
tain even more optimistic prediction for the first part of 
the section IV. 



APPENDIX B 

The long-range reconstructions of the sunspot number 
fluctuations (see, for instance, Refs. [3,0]) allow us to 
look on the solar transitional dynamics from a more gen- 
eral point of view. In figure 18 we show a spectrum of 
such reconstruction for the last 11,000 years (the data, 
used for computation of the spectrum, is available at 
[3]). The spectrum was computed using the autore- 
gressive model (with maximum entropy method [sH]). A 
semi-logarithmical representation was used in the figure 
to show an exponential law 



E{f) 



m 



The straight line is drawn in Fig. 18 to indicate the expo- 
nential law Eq. (Bl). Slope of the straight line provides 
us with the characteristic time scale = l//e — 176±7y. 
The exponential decay of the spectrum excludes the pos- 
sibility of random behavior and indicates the chaotic be- 
havior of the time series [4^,[3]. It is well known that 
low-order dynamic (deterministic) systems have as a rule 
exponential decay of E{f) (see, for instance, [il],[43]). 
As for infinite dimensional dynamic systems with chaotic 
attractors it is interesting to compare Fig. 18 with fig- 
ure 3 of the Ref. [48|. It should be noted that the 176y 
period is the third doubling of the period 22y. The 22y 
period corresponds to the Sun's magnetic poles polarity 
switching (see also Appendix C). 

The exponential spectrum can be also produced by a 
series of Lorentzian pulses with the average width of the 
individual pulses equals to r (though, the distribution of 
widths of the pulses should be fairly narrow to result in 
the exponential spectrum). 

In Fig. 18 a local maximum corresponding to the fre- 
quency fe and its first harmonics have been indicated 
by arrows. It should be noted that the harmonic 2/e 
corresponds to the well known period T ~ 88y (see, for 
instance, Ref. [i^), and the harmonic 3/e corresponds 
to the period T ~ 60y (cf. section "A FORECAST"). 

APPENDIX C 

It is also useful to consider the above results in content 
of multicentennial and millennial timescales. Figure 19 
shows a reconstruction of Northern Hemisphere tempera- 
tures for the past 2,000 years (the data for this figure were 



8 

6 

<C^ 4 

UJ 

C 2 



E~exp(-f/g 




0.1 



0.2 0.3 



0.4 



f my'] 



FIG. 18: Spectrum of the sunspots number fluctuations in 
the In- linear scales (the reconstructed data for the last 11,000 
years have been taken from [ij). The straight line is drawn 
to indicate the exponential law Eq. (Bl). 
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FIG. 19: A reconstruction of Northern Hemisphere temper- 
ature ano maly for the past 2,000 years (the data was taken 
from Ref. 



taken from Ref. [50|). This multi-proxy reconstruction 
was performed by the authors of Ref. [5l| using combina- 
tion of low-resolution proxies (lake and ocean sediments) 
with comparatively high-resolution tree-ring data. Fig- 
ure 20 shows a power spectrum of the data set calculated 
using the maximum entropy method (as for the spectra 
shown in Figs. 16 and 18) in the frames of the autore- 
gressive model, because it provides an optimal spectral 
resolution even for small data sets (see also [47|). The 
spectrum exhibits a rather wide peak indicating a peri- 
odic component with a period around 22 y, and a broad- 
band part with exponential decay. A semilogarithmical 
plot was used in Fig. 20 (cf Fig. 18) in order to show 
the exponential decay more clearly (at this plot the ex- 
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FIG. 20: Spectrum of the data, shown in Fig. 19, in semilog- 
arithmical scales. The straight line indicates the exponential 
decay Eq. (Bl). 
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FIG. 22: Spectrum of the data, shown in Fig. 19, in In-In 
scales. The dashed straight line indicates the Kolmogorov- 
like spectrum: E{f) ~ f~^^^ ■ The high-frequency part has 
been cutted in order to show the low- frequency part. 




FIG. 21; Spectrum of the chaotic fluctuations of the x- com- 
ponent for the Rossler system (a = 0.15, h — 0.20, c — 10.0.) 



ponential decay corresponds to a straight line). Both 
stochastic and deterministic processes can result in the 
broad-band part of the spectrum, but the decay in the 
spectral power is different for the two cases. The ex- 
ponential decay indicates that the broad-band spectrum 
for these data arises from a deterministic rather than a 
stochastic process. Indeed, for a wide class of determin- 
istic systems a broad-band spectrum with exponential 
decay is a generic feature of their chaotic solutions (see. 



for instance, Refs. [45| , [47| , [48| and Appendix B) 



In order to illustrate this and another significant fea- 
ture of the chaotic power spectra we show in figure 21 a 
power spectrum for the Rossler system [Hij 

dx , , dy dz ,^^\ 

-^ = -(y + ^); "^=^ + "y; —^b + xz-cz [CI) 

chaotic solution, where a, b and c are parameters. In 



this figure one can see a typical picture: a narrow-band 
peak (corresponding to the fundamental frequency of the 
system) in a low-frequency part and a broad-band expo- 
nential decay in a high-frequency part of the spectrum. 
Nature of the exponential decay of the power spectra 
of the chaotic systems is still an unsolved mathemat- 
ical problem. A progress in solution of this problem 
has been achieved by the use of the analytical contin- 
uation of the equations in the complex domain (see, for 
instance, [53|,[5J|)- ^^^^ approach the exponential de- 
cay of chaotic spectrum is related to a singularity in the 
plane of complex time, which lies nearest to the real axis. 
Distance between this singularity and the real axis deter- 
mines the rate of the exponential decay. If parameters of 
the dynamical system periodically fluctuate around their 
mean values, then at certain (critical) intensity of these 
fluctuations an additional singularity (nearest to the real 
time axis) can appear. Distance between this singularity 
and the real axis is determined by period of the peri- 
odic fluctuation of the system's parameters. Therefore, 
exponential decay rate of the broad-band part of the sys- 
tem spectrum equals the period of the parametric forcing. 
The chaotic spectrum provides two different characteris- 
tic time-scales for the system: a period corresponding to 
fundamental frequency of the system, Tjun , and a period 
corresponding to the exponential decay rate, Te = l//e 
(cf Eq. (Bl)). The fundamental period T/„„ can be esti- 
mated using position of the low-frequency peak, while the 
exponential decay rate period Tg = 1 //g can be estimated 
using the slope of the straight line of the broad-band part 
of the spectrum in the semilogarithmical representation 
(Figs. 20 and 21). From Fig. 20 we obtain T/„„ ~ 22±2y 
and Te ~ 11 ± ly (the estimated errors are statistical 
ones). Thus, the solar activity (SSN) period of 11 years 
is still a dominating factor in the chaotic temperature 
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fluctuations at the millennial time scales , although it is 
hidden for linear interpretation of the power spectrum 
(cf Ref. [5^). In the nonlinear interpretation the addi- 
tional period T^„„ ~ 22y might correspond to the funda- 
mental frequency of the underlying nonlinear dynamical 
system. It is surprising that this period is close to the 
22y period of the Sun's magnetic poles polarity switching 
(cf A ppe ndix B). It should be noted that the authors of 
Rcf. [56( found a persistent 22y cyclicity in sunspot ac- 
tivity, presumably related to interaction between the 22y 
period of magnetic poles polarity switching and a relic 
solar (dipole) magnetic field. Therefore, one cannot rule 
out a possibility that the broad peak, in a vicinity of fre- 
quency corresponding to the 22y period, is a quasi-/mear 
response of the global temperature to the weak periodic 
modulation by the 22y cyclicity in sunspot activity. I.e. 
strong enough periodic forcing results in the non-linear 
(chaotic) response whereas a weak periodic forcing results 
in a quasi-periodic response. 

Finally, figure 22 shows the same power spectrum as in 
Fig. 20 but in In- In scales, where the dashed straight line 
indicates the Kolmogorov-like behavior: E{f) ^ /^^^'^, 
in the lower-frequency part of the spectrum (see section 
'A Forecast' and Fig. 16 for relation between this type 
of spectrum and galactic turbulence). 
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